home *** CD-ROM | disk | FTP | other *** search
/ The 640 MEG Shareware Studio 2 / The 640 Meg Shareware Studio CD-ROM Volume II (Data Express)(1993).ISO / pascal / gaussj.zip / GAUSS.PAS < prev    next >
Pascal/Delphi Source File  |  1985-04-03  |  3KB  |  169 lines
  1. program gaus;    { -> 75 }
  2. { pascal program to perform simultaneous solution by Gaussian elimination }
  3. { procedure GAUSS is included }
  4.  
  5. const    maxr    = 8;
  6.     maxc    = 8;
  7.  
  8. type    ary    = array[1..maxr] of real;
  9.     arys    = array[1..maxc] of real;
  10.     ary2s    = array[1..maxr,1..maxc] of real;
  11.  
  12. var    y    : arys;
  13.     coef    : arys;
  14.     a    : ary2s;
  15.     n,m    : integer;
  16.     first,
  17.     error    : boolean;
  18.  
  19. external procedure cls;
  20.  
  21. procedure get_data(var a: ary2s;
  22.            var y: arys;
  23.          var n,m: integer);
  24.  
  25. { get values for n and arrays a,y }
  26.  
  27. var    i,j    : integer;
  28.  
  29. begin
  30.   writeln;
  31.   repeat
  32.     write('How many equations? ');
  33.     readln(n);
  34.     if not first then cls else first:=false;
  35.     m:=n
  36.   until n<maxr;
  37.   if n>1 then
  38.     begin
  39.       for i:=1 to n do
  40.     begin
  41.       writeln('Equation',i:3);
  42.       for j:=1 to n do
  43.         begin
  44.           write(j:3,':');
  45.           read(a[i,j])
  46.         end;
  47.       write(',C:');
  48.       read(y[i]);
  49.       readln    { clear line }
  50.     end;
  51.       writeln;
  52.       for i:=1 to n do
  53.     begin
  54.       for j:=1 to m do
  55.         write(a[i,j]:7:4);
  56.       writeln(':',y[i]:7:4)
  57.     end;
  58.      writeln
  59.     end        { if n>1 }
  60. end;    { procedure get_data}
  61.  
  62. procedure write_data;
  63.     { print out the answeres }
  64.  
  65. var    i    : integer;
  66.  
  67. begin
  68.   for i:=1 to m do
  69.     write(coef[i]:9:5);
  70.   writeln
  71. end;        { write_data }
  72.  
  73. procedure gauss
  74.            (a    : ary2s;
  75.         y    : arys;
  76.         var coef    : arys;
  77.         ncol    : integer;
  78.         var error    : boolean);
  79.  
  80. { matrix solution by Gaussian Elimination }
  81.  
  82. var
  83.     b    : ary2s;    { work array, nrow,ncol }
  84.     w    : arys;        { work array, ncol long }
  85.     i,j,i1,k,
  86.     l,n    : integer;
  87.     hold,sum,
  88.     t,ab,big: real;
  89.  
  90. begin
  91.   error:=false;
  92.   n:=ncol;
  93.   for i:=1 to n do
  94.     begin    { copy to work arrays }
  95.       for j:=1 to n do
  96.     b[i,j]:=a[i,j];
  97.       w[i]:=y[i]
  98.     end;
  99.   for i:=1 to n-1 do
  100.     begin
  101.       big:=abs(b[i,i]);
  102.       l:=i;
  103.       i1:=i+1;
  104.       for j:=i1 to n do
  105.     begin        { search for largest element }
  106.       ab:=abs(b[j,i]);
  107.       if ab>big then
  108.         begin
  109.           big:=ab;
  110.           l:=j
  111.         end
  112.          end;
  113.     if big=0.0 then error:= true
  114.     else
  115.       begin
  116.     if l<>i then
  117.       begin
  118.         { interchange rows to put largest element on diagonal }
  119.         for j:=1 to n do
  120.           begin
  121.         hold:=b[l,j];
  122.         b[l,j]:=b[i,j];
  123.         b[i,j]:=hold
  124.           end;
  125.           hold:=w[l];
  126.           w[l]:=w[i];
  127.           w[i]:=hold
  128.         end;    { if l<>i }
  129.       for j:=i1 to n do
  130.         begin
  131.           t:=b[j,i]/b[i,i];
  132.           for k:=i1 to n do
  133.         b[j,k]:=b[j,k]-t*b[i,k];
  134.           w[j]:=w[j]-t*w[i]
  135.         end    { j-loop }
  136.       end    { if big }
  137.     end;    { i-loop }
  138.       if b[n,n]=0.0 then error:=true
  139.       else
  140.     begin
  141.       coef[n]:=w[n]/b[n,n];
  142.       i:=n-1;
  143.       { back substitution }
  144.       repeat
  145.         sum:=0.0;
  146.         for j:=i+1 to n do
  147.           sum:=sum+b[i,j]*coef[j];
  148.         coef[i]:=(w[i]-sum)/b[i,i];
  149.         i:=i-1
  150.       until i=0
  151.     end;    { if b[n,n]=0 }
  152.       if error then writeln(chr(7),'ERROR: Matrix is singular')
  153. end;    { GAUSS }
  154.  
  155. begin        { MAIN }
  156.   first:=true;
  157.   cls;
  158.   writeln;
  159.   writeln('Simultaneous solution by Gauss elimination');
  160.   repeat
  161.     get_data(a,y,n,m);
  162.     if n>1 then
  163.       begin
  164.     gauss(a,y,coef,n,error);
  165.     if not error then write_data
  166.       end
  167.   until n<2
  168. end.
  169.